%reading the data file
close all
clear all
clc
DataFile = csvread('sc4002e0-rec_data.txt',1,0);
%DataFile
Time = DataFile(:,1);
EEG1 = DataFile(:,2);
%EEG1
EEG2 = DataFile(:,3);

Interval = 0.01;

%no of epochs
EpochLength = int32(30/Interval);%idivide(30,Interval,'floor');               
n = int32(length(Time)/EpochLength) - 1;%idivide(length(Time),EpochLength,'floor');
%n
data = EEG1;
%EpochLength
class = zeros(1, n);
for i = 1 : n
    
    startEpoch = (i-1)*EpochLength + 1;
    endEpoch = i*EpochLength;
    epoch = data(startEpoch : endEpoch);
    
    [Freq Amp] = FreqAmplitudeGraphInEpoch(epoch);
    
    avgAmp = Average(epoch);
    %avgAmp
    if  avgAmp > 10 && avgAmp < 30
        class(i) = 1; %awake
        continue;
    end
    
    %calculate avegrage amplitude in 2-7 Hz
    avgAmp2to7 = Average(BandPassFilter(epoch, 2, 7, Interval));
    
    if avgAmp2to7 > avgAmp && avgAmp2to7 > 100 && avgAmp2to7 < 200
        class(i) = 2; %stage1 sleep;
        continue;
    end
    
    %find signals with frequency in range 0-2 Hz
    freq2to7 = BandPassFilter(epoch, 0.000001, 2, Interval);
    
    %Now to calculate the percentage of signals having amp > 75microV
    %within epoch we split the epoch in intervals of 0.5s and calculate avg
    %amplitude of over that interval
    
    count = 0;
    for j = 1 : length(freq2to7)      
        if freq2to7(j) > 75
            count = count + 1;
        end 
    end
    if (count/EpochLength >= .2 && count/EpochLength <= .5)
        class(i) = 4;   %stage3
        continue;
    elseif count/EpochLength >= .5 
        class(i) = 5;   %stage4
        continue;
    end 
end
%class
count = 1;
%length(data)
for i = 1 : n
    
    startEpoch = (i-1)*EpochLength + 1;
    endEpoch = i*EpochLength;
    epoch = data(startEpoch : endEpoch);
    %length(epoch)
    %i
    %startEpoch
    %endEpoch
    for j = 1 : 0.5/Interval : EpochLength 
        endIndex = min(endEpoch, j + 0.5/Interval - 1);
        %j
        %endIndex
        %find signals with frequency in range 12-14 Hz
        freq12to14 = BandPassFilter(epoch(j : endIndex), 12, 14, Interval);
        if ~isempty(freq12to14)
            Spindle(count) = 1;
        else
            Spindle(count) = 0;
        end
        count = count + 1;
    end
end
%class
%Identifying sleep spindles and classifying as stage 1 or 2 based on that

for i = 1 : length(Spindle)
    j = i;
    if Spindle(i) == 1
        for j = i : min(length(Spindle), i + 359)
            if Spindle(j) == 1
                classTemp = 2;
                break;
            else
                classTemp = 1;
            end
            
        end
        if class(idivide(j, EpochLength, 'ceil')) == 0
            class(idivide(j, EpochLength, 'ceil')) = classTemp;
        end
    end
    i = j+1;
end
class
%calculating accuracy

%read hypnogram
hyp = csvread('sc4002e0-hyp_data.txt', 1, 0);
col1 = hyp(:,1);
col2 = hyp(:,2);
%col1(2)
%Convert to epoch based
count = 0;

for i = 1 : length(col1)
    %col1(i)
    epoch = int32(col1(i)/30);
    %epoch
    if epoch > length(class) - 1
        break;
    elseif class(epoch + 1) == int32(col2(i))
        count = count + 1;
    end
end

Accuracy = count/length(class)
